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Abstract 

We provide a status report on the calculation of the Equation of State (EoS) of QCD at finite 
temperature using lattice QCD. Most of the discussion will focus on comparison of recent results 
obtained by the HotQCD and Wuppertal-Budapest (W-B) collaborations. We will show that very 
significant progress has been made towards obtaining high precision results over the temperature 
range of T — 150 — 700 MeV. The various sources of systematic uncertainties will be discussed 
and the differences between the two calculations highlighted. Our final conclusion is that the 
lattice results of EoS are getting precise enough to justify being used in the phenomenological 
analysis of heavy ion experiments at RHIC and LHC. 
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1. The Road to Precision Lattice QCD Calculations 

One of the goals of simulations of lattice QCD is to provide a precise non-perturbative 
determination of the EoS of QCD over the temperature range 150 — 700 MeV that is 
being probed in experiments at RHIC and the LHC. The EoS, along with the transition 
temperature T c and transport coefficients such as shear viscosity, are crucial inputs into 
phenomenological hydrodynamical models used to describe the evolution of the quark 
gluon plasma (QGP). In this talk I will mainly review the two recent and most complete 
calculations of the EoS by the HotQCD [1] and Wuppertal-Budapest [2] collaborations. 

Simulations of lattice QCD at finite temperature are carried out on a 4-D hypercube 
of size aN T x aNs x aNs x aNs where a is the lattice spacing usually denoted in units 
of GeV -1 or fermi. The spatial size aNs is taken large enough so that finite volume 
corrections are under control and small. Past calculations show that for finite temperature 
simulations with Ns/N T — 4 these corrections are smaller than statistical errors for 
T<3T C . For higher temperatures larger Ns may be required. In QCD, the gauge coupling 
ft = 6/g 2 is related to a by dimensional transmutation and the continuum theory is 
recovered in the limit a->0or equivalently g — > or j3 — > oo. To provide a perspective 
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on how fine the current lattice simulations are note that N T — 10 corresponds to a « 0.1 
fcrmi or ~ 2 GeV^ 1 at the transition temperature T ~ 200 MeV. 

The second set of control parameters (inputs) in the simulations are the quark masses. 
The results discussed in this review are for 2 + 1 flavors, i.e., two flavors of degenerate 
u and d quarks and a heavier strange quark. In nature 2m s /(m u + m<j) w 27.5 and 
to = (m u + m<j)/2 ~ 3.5 MeV and m s <~ 90 MeV. The small value of to provides the 
key computational challenge because the most time consuming part of the simulations 
is the inversion of the Dirac operator, a very large sparse matrix. This inversion is done 
using interative Krilov solvers that have critical slowing down in the limit to — > 0. The 
computational cost increases as to~ 7 / 2 or faster and calculations become very expensive 
with decreasing quark mass. For this reason 2 + 1 flavor simulations are typically done 
fixing the strange quark mass to its physical value and simulating at a number of values 
of m/m s from which extrapolations to the physical to are made. Recent calculations 
by the W-B and HotQCD collaborations show that the computer power has reached a 
stage where simulations can be done close to, or directly at, physical to. Thus, in the 
state-of-the-art simulations, this source of systematic error (extrapolation in to) is now 
under good control. 

The tuning of the set of parameters {g,m s ,m} is done as follows. One first fixes the 
ratio m s /m, ideally to m s /m = 27.5. Then for a judiciously chosen value of g, zero 
temperature simulations are done to measure two independent physical quantities whose 
values are experimentally measured or well determined, and one of which is sensitive to 
strange quark mass, for example the K-meson mass Mk and the pion decay constant f v 
(or fx and M n ). The value of m s is tuned until the lattice results for these two quantities 
match their physical values. This fixes a and m s . Now depending on how finely one wants 
to scan in T (or a) a new value of g is chosen and the value of m s is again tuned to 
reproduce the observables, thus determining the new a keeping to s /to fixed. This process 
generates a set of {g, m s ,m} values for which, by construction, the physics (defined by 
matching lattice Mk and f v to physical values) is fixed. This line in the {g, m s ,m} space, 
since m s is tuned to the physical value and m s /m is fixed, is called a line of constant 
physics (LoCP) . The utility of simulating along LoCP is to reduce the three dimensional 
space of input parameters to a line along which only the lattice spacing is changing. 
This procedure provides better control over taking the continuum limit. The extent to 
which {g,m s } would have varied had one chosen two different physical quantities, say 
Mg S and Mjy, is a measure of variations in discretization errors along different LoCP. The 
emphasis of the lattice community is to present results extrapolated to a — > in order to 
remove these discretization errors, which are proportional to powers of a, by simulating 
at a number of values of a. In case of finite temperature calculations we extrapolate 
results at fixed T to iV r — > oo by simulating at a number of values of N T . 

The lattice size in the Euclidean time direction defines the temperature of the system by 
the relation T = l/aN T . The scale a for fixed {g, m s ,fn} is the same for zero-temperature 
and finite temperature lattices. Thus, knowing a corresponding to a given g and quark 
masses uniquely determines T for a given N T . An important consequence of the fact that 
g or equivalently a or T is the single parameter that controls lattice simulations is that 
only one thermodynamic quantity can be determined, which for the extraction of EoS is 
the trace anomaly I/T 4 = (e — 3p)/T 4 , where I is called the integration measure. 

The above approach for scanning in T is called the fixed N T approach. A second 
approach that I will not cover in detail and which is being pursued by the WHOT 
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collaboration [3] with improved Wilson fermions is called the fixed a approach. In this 
approach, for a given a (the same process is used for fixing {<?, m s , m}) one simulates on 
a number of different N T lattices to scan in T . The advantage of this approach is that 
only a single zero-temperature matching calculation, needed to carry out subtractions of 
lattice artifacts in finite T data, is required for each a. The weakness of this approach 
is that the scan in T is limited by the coarseness and range of N T values possible, i.e. 
N T = 6, 8, 10, 12, 14, before one runs out of computational power. The recent results using 
the fixed a approach by the WHOT collaboration [3] are very encouraging and I refer 
interested readers to their paper for details. 

2. Taste Symmetry breaking with Staggerd fermions 

In the naive discretization of the Dirac action there automatically is a 2 4 doubling 
of flavors. In the staggered fcrmion approach, using a lattice symmetry called "spin 
diagonalization" , the degeneracy is reduced from 16 to 4 by placing a single degree of 
freedom at each lattice site. Under this construction a 2 4 hypercube is the basic unit cell 
that reduces to a point in the continuum limit. The 16 degrees of freedom in the cell 
represent, in the continuum limit, four identical copies (called taste) of Dirac fermions. On 
the lattice this four-fold degeneracy gives rise to a proliferation of particles propagating 
in the QCD vacuum, for example, there are 16 pions distinguished by their taste rather 
than just one. If taste symmetry was unbroken at finite a, the 4-fold degeneracy could 
be handled by just dividing results by the appropriate degeneracy factor. Problems arise 
because this degeneracy is broken at finite a and one does not know a priori how large this 
systematic error is and how it effects simulations of QCD thermodynamics in particular. 
The most common approach to quantify this effect is to study the difference in the 
masses of the 16 pions, and determine how these differences vanish when lattice results 
are extrapolated to the continuum limit. 

The HotQCD collaboration has studied the consequences of taste symmetry breaking 
utilizing three versions of improved staggered fermions — the asqtad, p4 and HISQ/tree 
formulations [4]. In Fig. 1 we show preliminary HotQCD HISQ/tree results for M% — Mq 
where M G is the Goldstone pion mass versus a, and also compare with results from the 
stout and asqtad actions. The large spread in masses that increases on coarser lattices 
shows that taste symmetry is indeed badly broken [5] . Based on such studies the conclu- 
sion is that at any given a, the taste breaking is least in HISQ/tree followed by stout, 
asqtad, and p4 actions. One consequence of this discretization error for thermodynamics 
is that the contribution of any state, for example the pion, is not just from the lowest 
taste state (the Goldstone pion) but is some weighted average of the 16 pions (or the 
appropriate multiplets for other states). Thus, at any given a the effective masses of all 
hadrons arc larger than the desired ground state value. The magnitude of the effect is 
shown in the right panel of Fig. 1 which plots the root mean squared mass of the 16 pion 
states corresponding to a Goldstone pion mass of 140 MeV. Results at low temperatures, 
T < 150 MeV, and on small N T lattices are most susceptible to this discretization error. 
Taste breaking also puts a caveat on the above described tuning of quark masses: setting 
m s /m = 27.5 does not guarantee that simulations have been done at the physical light 
quark masses. Simulating at a number of m s /m values to control the chiral extrapolation 
followed by taking the continuum limit using a number of N T lattices provides the best 
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Fig. 1. The splitting between the 16 pions due to taste symmetry breaking for HISQ/tree, stout and 
asqtad actions (left panel). The right panel shows the root mean squared mass corresponding to a 
Goldstone pion mass of 140 MeV. Note that the HISQ action used so far for thermodynamics does not 
include tadpole improvement of the gauge action and is therefore called HISQ/tree in [4]. 



understanding of systematic errors and for obtaining physical results. 

A second issue with staggered simulations that include the strange quark (or in future 
the charm quark) is the need to take the fourth root of the determinant to compensate 
for the 4-fold degeneracy. Creutz [6] claims that this "rooting" is a fundamental flaw 
of the staggered formulation, however, its effect in current simulations may be small 
since the strange quark mass is large, whereas the review by Sharpe [7] (covering a large 
body of work) shows that while the staggered formulation may be ugly, it gives physical 
results in the continuum limit. (For degenerate u and d quarks, this rooting problem 
is overcome because of a lattice symmetry whereby the square root of the determinant 
can be written as the determinant on even (or odd) sites.) In a perfect world one would 
like to use a Wilson-like action that maintains the continuum flavor structure and chiral 
symmetry, i.e., domain wall or overlap fermions, however, to date most thermodynamical 
simulations use staggered fermions for two reason: they are much faster (10 — 20x) to 
simulate than even simple Wilson fermions and because of a residual chiral symmetry 
that protects the Goldstone pion. 



3. The Trace Anomaly 

The results from the HotQCD [1] [4] and W-B [2] collaborations for I/T 4 = (e-3p)/T 4 , 
the single thermodynamic quantity calculated on the lattice, are shown in Fig. 2. 

Before making detailed comparisons it is important to stress here, and applicable to all 
discussion that follows, that the HotQCD results do not yet incorporate extrapolation 
to the physical quark mass or the continuum limit. The most extensive data are for 
N r = 6, 8 and m s /fn = 10 with new ongoing calculations at m s /fn = 20 and N T = 12. 
The W-B results at N T = 6, 8 have been obtained for number of values of m s /m including 
at m s /m = 28.15 where the final physical value is quoted. (Recall, however, the caveat 
about the uncertainty in locating the physical value of m due to taste breaking.) The 
W-B data at N T = 10 and 12 are more limited in T values and are at m s /m — 28.15 
only. W-B define their N T — 8 results to represent the continuum value. Data with 
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Fig. 2. Integration measure I/T 4 from W-B (left) and preliminary HotQCD (right) comparing different 
actions. The curve (right panel) is a a paramctrization of e — 3p/T 4 that is based on HRG and the 
p4/asqtad lattice data at high temperatures [8]; and the W-B data are shown in purple for comparison. 



N T = 10 at T < 365 MeV (red points in left panel of Fig. 2) and the three points (green 
points in Fig. 2) with N T = 12 provide a consistency check since they show no significant 
discretization effects relative to N T = 8 data. 

The overall form of the results by the two collaborations is similar. There are, however, 
two significant differences between the HotQCD and W-B data. The first is the value of 
I at the peak, ~ 5.8 versus 4.1, and the peak in the W-B data is shifted to lower T by 
about 20 MeV. W-B collaboration attribute these differences to the lack of extrapolation 
of HotQCD data in quark mass and a, i.e., residual discretization errors. Preliminary 
HotQCD results with N T — 8 lattices using HISQ/tree fermions (also shown in Fig. 2) 
give a similar peak height and position as the asqtad action (the p4 results are higher 
but decreasing with N T ). The agreement between asqtad and HISQ/tree and since the 
HISQ/tree action is the more improved than stout, has smaller discretization errors and 
less taste breaking, it is not clear if, today, we have a simple resolution of the difference. 
The forthcoming results with HISQ/tree action on N T = 8 and 12 lattices being simulated 
by the HotQCD collaboration should help clarify these issues. 



4. Pressure, Energy Density, Entropy and Speed of Sound 



The pressure p can be determined from the trace anomaly using the following relations: 



I e^(T) e - 3p ( P ^ 
P(T) p(T ) } 



T 4 T 4 ~ J dt t 5 ( 2 ) 

To 

The results for pressure and energy density are summarized in Figs. 3 and 4. To obtain 
pressure p there are two issues that need to be addressed when carrying out the inte- 
gration in Eq. 2. The first is to construct a smooth function that represents the lattice 
data for (e — 3p)/T 4 over the whole range of T since I/T 4 has been calculated only at a 
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Fig. 3. HotQCD collaboration [1] results for e/T 4 and 3p/T 4 on N T = 8 lattices with m s /m = 10. 

finite number of values of T. The second is the choice of To above which 9 mm (To)/Tq is 
well-determined and at which point p(T ) can be estimated reliably. 

The HotQCD collaboration has investigated a number of ansatz for parameterizing 
e w (T )/T 4 and find that the results for p do not vary significantly. The uncertainty due 
to the ansatz is shown by the error bars on p/T A at T = 275, 540 MeV in Fig. 3. The 
W-B collaboration uses a variant of the method — they parameterize the pressure itself 
and then evaluate its derivatives to match to /. The band in Fig 4 show the uncertainty. 

The second issue is more significant. The systematic errors in lattice data grow as T is 
lowered and are expected to be large below T = 150 MeV. At the same time p(T = 150) 
is not negligible and o priori unknown. One approach is to use the hadron resonance gas 
(HRG) model for p(T = 150 MeV). This requires that there be reasonable agreement 
between the HRG and lattice values at T = 150 MeV. The HotQCD data approaches the 
HGR from below and at T = 150 MeV there is a significant difference. Another approach 
is to use To = 100 MeV where there is more confident in the HRG value but then one 
has to confront the uncertainty in matching and parameterizing mm (To)/Tq between 
T = 100 — 150 MeV. The HotQCD collaboration use p = at T = 100 MeV for their 
central value and the HRG value to estimate the uncertainty whose magnitude is shown 
by the black square on e/T 4 at T = 550 MeV in Fig. 3. The W-B collaboration show 
that a modified "lattice" HRG calculation, taking into account taste breaking in pion 
and kaon states, fits the lattice data between T = 100 — 140 MeV. Nevertheless, they 
choose p(T = lOOAfeV, m s = m) = for the normalization. 

Once / and p are determined the energy density is given by e/T 4 = 7/T 4 + 3p/T 4 , 
entropy by s = (e + p) /T and the speed of sound c s by 

c 2 = ^P = c d (p/£) P 
s de de e 

A comparison of results for e, p and s is shown in Figures 3, 4, and 5. 

The W-B collaboration apply two corrections to the estimate for p(T). The first is to 
guarantee that the lattice results for each N T match the continuum Stefan-Boltzmann 
value at T = oo. To do this they construct the ratio of the continuum Stefan-Boltzmann 
value for p to its free-field (T = oo) lattice value (the continuum integrals are replaced by 
lattice sums appropriate for each N T ) and then correct the lattice data at all T by this 
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Fig. 4. Results for e/T 4 and p/T 4 from the W-B collaboration [2]. To compare results for p note that 
the W-B data are for p/T 4 whereas the HotQCD results in Fig 3 are for 3p/T 4 . 
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Fig. 5. Comparison of entropy density obtained by W-B (left) and HotQCD (right) collaborations. The 
HotQCD data are for m s /fn = 10. 

ratio. This ratio is large, 1.517 and 1.283 for N T — 6 and 8 respectively. Since N T = 8 
data are used to define the continuum estimate, this correction is too large to justify 
on the basis of a tree-level improvement of lattice observables (lattice operators used to 
probe the physics). Furthermore, this correction is also applied to /, e and s. The second 
correction made by the W-B collaboration is to shift upward their results for p/T 4 by 
half the difference, 0.06, between the lattice and HRG estimates at T — 100 MeV. Again, 
to me, this is not a well-motivated correction of data. Future simulations and better 
understanding of the low T region will hopefully alleviate the need for such corrections. 

A comparison of results for the speed of sound are shown in Fig. 6. The fundamental 
quantity needed to calculate it is p/e as shown in Eq. 3. Two features in the data are 
worth commenting on. First, data in Fig. 4 show that, in the transition region from 
hadronic matter to QGP, the energy density is changing more rapidly than the pressure. 
This implies that c s should show a dip in the transition region as is indeed observed. 
Second, c s rises quickly after the transition region and reaches close to the relativistic 
Boltzmann gas value of 1/3 by T ~ 400 MeV. 

5. Prospects for improvement in the EoS of QCD in the near future 

Significant progress has been made in determining the EoS using lattice QCD in the 
last three years. The current lattice results for the EoS have already given the Heavy Ion 
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Fig. 6. Comparison of the speed of sound c^. The W-B data (left panel) are plotted versus T. The 
HotQCD data for p/e and versus the energy density e using the fits to e and p (right panel). 

community a much better understanding of the dynamics of the QGP. Lattice estimates 
are now being used in hydrodynamic analysis of the evolution of the QGP in experiments 
at the RHIC at BNL and at the LHC. 

There are a number of ways in which both HotQCD and W-B collaborations are 
improving their estimates: 

- The HotQCD collaboration will present their estimates of the continuum values with 
HISQ/tree action on N T = 6, 8 and 12 lattices and with m s /m = 10 and 20. 

- Both collaborations will include the charm quark to provide results with (2 + 1 + 1) 
dynamical flavors. Preliminary partially-quenched estimates suggest that the charm 
quark contribution starts to become large at above T ~ 300 MeV. 

- To fully control finite volume effects at T > 500, both collaborations will simulate on 
larger spatial volumes, larger than N S /N T = 4. 

- In addition to simulations with staggered fermions, simulations with improved Wilson 
and domain wall fermions are maturing [3] [9] and will provide independent checks of 
the staggered results in, hopefully, the near future. 

With these improvements a number of unresolved issues such as the location and height 
of the peak in I, control over systematic errors at low temperatures, and the impact of 
charm quark at high temperatures, should be addressed over the next couple of years. 
So stay tuned. 
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